1. Field of the Invention
The present invention is directed to three-dimensional surface profile imaging, and more particularly to a method and apparatus for three-dimensional imaging that uses color ranging to conduct surface profile measurement.
2. Description of the Related Art
A three dimensional surface profile imaging method and apparatus described in U.S. Pat. No. 5,675,407 (xe2x80x9cthe ""407 patentxe2x80x9d), the disclosure of which is incorporated herein by reference in its entirety, conducts imaging by projecting light through an optical filter, such as linear variable wavelength filter (LVWF), thereby projecting light having a known, spatially distributed wavelength spectrum on the objects being imaged. The LVWF is a rectangular optical glass plate coated with a color-filtering film that gradually varies in color, (i.e., wavelength). If the color spectrum of a LVWF is within the visible light region, one edge of the filter rectangle may correspond to the shortest visible wavelength (i.e. blue or violet) while the opposite edge may correspond to the longest visible wavelength, (i.e. red). The wavelength of light passing through the coated color-filtering layer is linearly proportional to the distance between the position on the filter glass where the light passes and the blue or red edge. Consequently, the color of the light is directly related to the angle xcex8, shown in FIG. 1, at which the light leaves the rainbow projector and LVWF.
Referring to FIG. 1 in more detail, the imaging method and apparatus is based on the triangulation principle and the relationship between a light projector 100 having the LVWF, a camera 102, and the object being imaged 104. As shown in FIG. 1, a triangle is uniquely defined by the angles theta (xcex8) and alpha (xcex1), and the length of the baseline (B). With known values for xcex8, xcex1, and B, the distance (i.e., the range R) between the camera 102 and a point Q on the object""s surface can be easily calculated. Because the baseline B is predetermined by the relative positions of the light projector 100 and the camera 102, and the value of xcex1 can be calculated from the camera""s geometry, the key to the triangulation method is to determine the projection angle, xcex8, from an image captured by the camera 102 and more particularly to determine all xcex8 angles corresponding to all the visible points on an object""s surface in order to obtain a full-frame 3D image in one snapshot.
FIG. 2 is a more detailed version of FIG. 1 and illustrates the manner in which all visible points on the object""s surface 104 is obtained via the triangulation method. As can be seen in the Figure, the light projector 100 generates a fan beam of light 200. The fan beam 200 is broad spectrum light (i.e., white light) which passes through the LVWF to illuminate one or more three-dimensional objects 104 in the scene with a pattern of light rays possessing a rainbow-like spectrum distribution. The fan beam of light 200 is composed of multiple vertical planes of light, or xe2x80x9clight sheetsxe2x80x9d, each plane having a given projection angle and wavelength. Because of the fixed geometric relationship among the light source 102, the lens of the camera 102, and the LVWF, there exists a one-to-one correspondence between the projection angle (xcex8) of the vertical plane of light and the wavelength (xcex) of the light ray. Note that although the wavelength variations are shown in FIG. 2 to occur from side to side across the object 104 being imaged, it will be understood by those skilled in the art that the variations in wavelength could also be made from top to bottom across the object 104 being imaged.
The light reflected from the surface of the object 104 is then detected by the camera 102. If a visible spectrum range LVWF (400-700 nm) is used, the color detected by the camera pixels is determined by the proportion of its primary color Red, Green, and Blue components (RGB). The color spectrum of each pixel has a one-tone correspondence with the projection angle (xcex8) of the plane of light due to the fixed geometry of the camera 102 lens and the LVWF characteristics. Therefore, the color of light received by the camera 102 can be used to determine the angle xcex8 at which that light left the rainbow light projector. Other spectrum ranges can also be used in similar fashion. The implementation is straightforward to those skilled in the art.
As described above, the angle a is determined by the physical relationship between the camera 102 and the coordinates of each pixel on the camera""s imaging plane. The baseline B between the camera""s 102 focal point and the center of the cylindrical lens of the light projector 100 is fixed and known. Given the value for angles xcex1 and xcex8, together with the known baseline length B, all necessary information is provided to easily determine the full frame of three-dimensional range values (x,y,z) for any and every visible spot on the surface of the objects seen by the camera 102.
As shown in FIG. 3, given the projection angle xcex8, the three-dimensional algorithm for determining the (x,y,z) coordinates of any surface spot Q(x,y,z) on a three-dimensional object is given below based on the following triangulation principle:                                           x            =                                          B                                                      f                    *                    ctg                    ⁢                                          xe2x80x83                                        ⁢                    θ                                    -                  u                                            *              u                                ,                      xe2x80x83                    ⁢                      y            =                                          B                                                      f                    *                    ctg                    ⁢                                          xe2x80x83                                        ⁢                    θ                                    -                  u                                            *              v                                ,                      xe2x80x83                    ⁢                      z            =                          xe2x80x83                        ⁢                                          B                                                      f                    *                    ctg                    ⁢                                          xe2x80x83                                        ⁢                    θ                                    -                  u                                            *              f                                      ⁢                  xe2x80x83                                    (        1        )            
As a result, the three dimensional imaging system described above can capture full-frame, high spatial resolution three-dimensional images using a standard camera, such as a charge coupled device camera, in real time without relying on any moving parts. Further, because the imaging system does not rely on a laser, it does not pose any hazard to the eyes when used in clinical applications. Also, because the wavelength of the light projected onto the object surface continuously varies, there is no theoretical limitation on the measurement accuracy that can be achieved by the system. The actual accuracy of a specific system will depend on system implementation and will be affected primarily by limiting factors such as the optical system design, the quality and resolution of the camera, light spectral emission of the light source projector; noise level and resolution of the frame grabber, calibration algorithms, and the three-dimensional imaging processing algorithms.
To avoid allowing the ambient light on the object being imaged from affecting the imaging results, the system may obtain an image of the object under normal light conditions before projecting the filtered light onto the object. The image obtained under normal light conditions is then subtracted from the image obtained under LVWF light conditions to eliminate the effects of the ambient light on the image.
Referring back to FIG. 1, the triangulation algorithm used in the imaging system is based on the following formula, with reference to FIG. 3:                               R          =                                                    sin                ⁢                                  xe2x80x83                                ⁢                θ                                            sin                ⁢                                  xe2x80x83                                ⁢                α                                      ⁢            B                          ,                            (        2        )            
where (xp, yp) is the location of the rainbow light projector, (xc, yc) is the location of imaging sensor, B is the baseline between the rainbow projector and the imaging sensor (CCD), xcex1=xcfx80xe2x88x92xcex8xe2x88x92xcex2, O is a surface point on the object in the scene, and R is the three-dimensional range, that is, the distance between (xc, yc) and O.
Note that all of the variables, xcex8, xcex1, and xcex2, in the equation 10 may introduce error in the three-dimensional range calculation. In the following error sensitivity analysis, considered with reference to FIG. 4, it is assumed that the coordinate of the camera""s focal point in the world coordinate system can be obtained precisely through camera calibration. The full derivative of the R is given by:                               ⅆ          R                =                                            [                                                                    sin                    ⁢                                          xe2x80x83                                        ⁢                    θ                    ⁢                                          xe2x80x83                                        ⁢                    cos                    ⁢                                          xe2x80x83                                        ⁢                    α                                                                              -                                              sin                        2                                                              ⁢                    α                                                  ⁢                b                            ]                        ⁢                          ⅆ                              xe2x80x83                            ⁢              α                                +                                    [                                                                    cos                    ⁢                                          xe2x80x83                                        ⁢                    θ                                                        sin                    ⁢                                          xe2x80x83                                        ⁢                    α                                                  ⁢                b                            ]                        ⁢                          ⅆ                              xe2x80x83                            ⁢              θ                                +                                    [                                                sin                  ⁢                                      xe2x80x83                                    ⁢                  θ                                                  sin                  ⁢                                      xe2x80x83                                    ⁢                  α                                            ]                        ⁢                          ⅆ              b                                                          (        3        )            
The dxcex1, dxcex8, and db are all functions of the position error of the light projector (dxp, dzp) and the color mismatch dx. Notice that:                               b          =                                                                      (                                                            x                      c                                        -                                          x                      s                                                        )                                2                            +                                                (                                                            z                      c                                        -                                          z                      s                                                        )                                2                                                    ,                  
                ⁢                  β          =                                                    -                                  tan                                      -                    1                                                              ⁢                                                                    z                    c                                    -                                      z                    p                                                                                        x                    c                                    -                                      x                    p                                                                        +                                          tan                                  -                  1                                            ⁢                                                z                  c                                                                      x                    c                                    -                  x                                                                    ,                  xe2x80x83                ⁢                  and          ⁢                      
                    =                                                    tan                                  -                  1                                            ⁢                                                                    z                    c                                    -                                      z                    p                                                                                        x                    c                                    -                                      x                    p                                                                        +                                          tan                                  -                  1                                            ⁢                                                z                  p                                                  x                  -                                      x                    p                                                                                                          (        5        )            
Therefore                               db          =                                                                      (                                                            x                                              p                        ⁢                                                  xe2x80x83                                                                                      -                                          x                      c                                                        )                                ⁢                                  dx                  p                                            +                                                (                                                            z                      p                                        -                                          z                      c                                                        )                                ⁢                                  dz                  p                                                      b                          ⁢                  
                ⁢                              d            ⁢                          xe2x80x83                        ⁢            β                    =                                                                      (                                                            x                                              cf                        ⁢                                                  xe2x80x83                                                                                      -                                          x                      p                                                        )                                ⁢                                  dz                  p                                            +                                                (                                                            z                      p                                        -                                          z                      c                                                        )                                ⁢                                  dx                  p                                                                    b              2                                      ⁢                  
                ⁢                              d            ⁢                          xe2x80x83                        ⁢            θ                    =                                                                                          (                                                                  x                                                  p                          ⁢                                                      xe2x80x83                                                                                              -                                              x                        c                                                              )                                    ⁢                                      dz                    p                                                  +                                                      (                                                                  z                        c                                            -                                              z                        p                                                              )                                    ⁢                                      dx                    p                                                                              b                2                                      -                                                            z                  p                                                                                            (                                              x                        -                                                  x                          p                                                                    )                                        2                                    +                                      z                    p                    2                                                              ⁢                              xe2x80x83                            ⁢              dx                                                          (        5        )            
To simplify the analysis, the effect of a mismatch error in the world coordinate, dx, caused by imaging sensor""s spectral noise can be indicated as follows:                               ⅆ          R                =                              [                                                            cos                  ⁢                                      xe2x80x83                                    ⁢                  θ                                                  sin                  ⁢                                      xe2x80x83                                    ⁢                  α                                            ⁢              b                        ]                    ⁢                                    z              p                                                                        (                                      x                    -                                          x                      p                                                        )                                2                            +                              z                p                2                                              ⁢                      xe2x80x83                    ⁢                      ⅆ            x                                              (        6        )            
The linear relationship between the mismatch error dx in the world coordinate and the mismatch error in CCD camera""s pixel coordinate, di, allows approximation of dx and di as if they are the same order of magnitude, resulting in the expression:
dx=k*dixe2x80x83xe2x80x83(7)
where k is a constant. From this expression, the range error in terms of the color mismatch error in CCD pixel coordinate can be expressed as follows:                               ⅆ          R                =                              [                                                            cos                  ⁢                                      xe2x80x83                                    ⁢                  θ                                                  sin                  ⁢                                      xe2x80x83                                    ⁢                  α                                            ⁢              b                        ]                    ⁢                                    z              p                                                                        (                                      x                    -                                          x                      p                                                        )                                2                            +                              z                p                2                                              *          k          *                      ⅆ            i                                              (        8        )            
In the imaging system shown in FIGS. 1 and 2, the continuously varying wavelength pattern obtained via the light projector 100 and the LVWF 101 is projected spatially onto surfaces of a three-dimensional object or scene 104. There is a one-to-one corresponding relationship between the wavelength of the projected light sheet, xcfx89, and its projection angle, xcex8. Based on this one-to-one relationship, a (xcfx89,xcex8) lookup table can be established for the color matching operation. An image sensor device (such as the CCD camera 102) that collects the reflection from object""s surfaces also provides wavelength information for each pixel in an image, enabling the determination of the projection angle xcex8 corresponding to each image pixel through the color matching operation based on the (xcfx89,xcex8) lookup table.
The accuracy of the color match operation, of course, affects the accuracy of the three-dimensional surface measurement. Referring to FIG. 5, the accuracy of the color match is significantly determined by the characteristics of the projected wavelength-varying illumination pattern. More particularly, the system has a signal to noise ratio (S/N=c), where the noise is a catch-all term that characterizes all undesirable effect in the system to reduce the fidelity of the spectral signal (i.e., wavelength). Due to the existence of the noise, mismatch error (denoted as dw) occurs.
As shown in FIG. 6, for the same magnitude level of dw, if we use a spatially varying projection pattern that has a greater derivative, a smaller mismatch error in the pixel coordinate, denoted as di, occurs due to the increased slope caused by the larger derivative:
di2 less than di1.xe2x80x83xe2x80x83(9)
Using Equ. (8) in Equ. (9), the expression becomes:
dR2 less than dR1.xe2x80x83xe2x80x83(10)
The above expressions suggest that mismatch errors can be reduced by using a spatially varying projection pattern that has a large variation rate and therefore a large derivative.
As noted above, the imaging system projects a light pattern that could be described as a xe2x80x9csingle rainbow projectionxe2x80x9d pattern because the pattern ranges over the visible light spectrum a single time. This type of pattern, however, has a variation rate that depends on numerous factors, including the field of view of the camera/image sensor, the distance between the light projector and the object being imaged, as well as other factors. As a result, it is difficult to adjust the spectral variation rate of the single rainbow projection pattern to fit the accuracy requirements of many three-dimensional imaging tasks.
There is a need for a method and system that can optimize the light pattern to provide higher measurement accuracy.
Accordingly, the present invention is directed to a system and method that generates a light pattern that contains multiple cycles of light spectra or sub-patterns. In one embodiment, the light pattern includes two or more cycles of spectra of visible or infrared light, arranged contiguously with each other. The light pattern may also be generated by sequentially emitting two or more sub-patterns.
In another embodiment of the invention, the multiple spectra light pattern is generated by projecting the light from the light source through a planar member, such as a slotted panel, that is movable to generate small phase shifts in the light pattern. The multiple patterns and/or sub-patterns used to illuminate the object being imaged improves resolution and measurement accuracy by reducing color mismatch error.